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A new method is presented for estimating frequency responses and their uncertainties 
from wind-tunnel data in real time. The method uses orthogonal phase-optimized multi- 
sine excitation inputs and a recursive Fourier transform with a least-squares estimator. The 
method was first demonstrated with an F-16 nonlinear flight simulation and results showed 
that accurate short period frequency responses were obtained within 10 seconds. The 
method was then applied to wind-tunnel data from a previous aeroelastic test of the Joined- 
Wing SensorCraft. Frequency responses describing bending strains from simultaneous 
control surface excitations were estimated in a time-efficient manner. 


Nomenclature 


Symbols 


a, b 

input and output amplitudes 

G{ju>) 

frequency response 

3 

imaginary number \/— 1 

K 

set of excitation frequency indices 

N 

number of data samples 

n f 

number of Fourier transform samples 

n u 

total number of excitation frequencies 

5ft 

real part 

T 

record length [s] 

t 

time [s] 

u 

input 

X 

regressor matrix 

y 

output 

A 

perturbation 

6 

control surface deflection [deg] 

e 

bending strain 

e 

unknown parameters 

s 

covariance matrix 


a 1 2 

variance 


input and output phase angles [rad] 

CO 

frequency [rad/s] 

INI 

magnitude 

L 

phase angle 

Superscripts 

T 

transpose 

-1 

inverse 

t 

complex conjugate transpose 
time derivative 


estimated value 
Fourier transform 

Subscripts 

AW, FW 

aft and forward wing 

L, R 

left and right side 

O, M, I 

outer, middle, and inner 
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I. Introduction 


T O meet increasing demands for fuel economy gain and noise pollution reduction through improved aero- 
dynamic efficiency, next-generation aircraft are expected to be light weight and highly sophisticated, as 
indicated by the Green Aircraft Challenge and the NASA Fundamental Aeronautics Program, for example. 
As a consequence of these concept aircraft incorporating lighter, larger, and more flexible components, aero- 
elastic interactions become significant and a substantial effort is needed to understand, model, and mitigate 
these effects. 

Frequency responses are often desired because they are non-parametric linear models that can be iden- 
tified from measured data and do not require prior knowledge of the model structure. In aeroelastic wind 
tunnel tests, each control surface on the test article is sequentially excited using a frequency sweep, which 
consists of a sinusoid waveform with a continuously increasing frequency. After the test has finished and 
data has been collected, frequency responses between each input and output are computed using Welch’s 
method or a similar technique; 1-3 however, this method has a number of drawbacks. Frequency sweeps 
continuously transition through the frequencies, which excites transient responses and limits the amount 
of excitation at any given frequency. The input amplitude is designed according to the output response 
amplitude at resonant frequencies, which results in lower signal-to-noise ratios when exciting other frequen- 
cies. For rigid-body models, long record lengths of 90-120 seconds of data are needed for good results. 2 In 
aeroservoelastic problems, there are numerous inputs and outputs due to the spatially distributed nature of 
structural dynamics, which can further lead to long test times. While several commercial packages exist to 
help automate the analysis, 2,4 numerous engineering judgments are needed to estimate frequency responses, 
including selection of window number, size, shape, and overlap, to mitigate problems due to finite record 
lengths, low excitations, and measurement noise. Finally, there is no direct measurement of uncertainty for 
the estimated frequency response; instead, the coherence is used, which evaluates the degree of linearity 
between an input and output. 

A new method for estimating frequency responses and their uncertainties efficiently and in real time is 
presented here. The method uses orthogonal phase-optimized multi-sine inputs to excite all control surfaces 
simultaneously, uniquely, and over a specified frequency range of interest. The Fourier transform converts 
output data into the frequency domain, where a least-squares technique extracts the relative magnitude and 
phase angle, as well as the associated uncertainties. A recursive version of the Fourier transform is also 
implemented to perform the estimation in real time as the data streams, rather than after the test has been 
conducted. This method removes many of the engineering judgments needed in frequency response estimation 
and directly estimates the uncertainty. Consequently, this method can change the way in which wind tunnel 
tests occur: data is streamed until frequency response estimates have reached prescribed uncertainty levels, 
at which time the test is stopped. This procedure would help diagnose problems and minimize time spent 
conducting wind tunnel tests. The method is applied here to previous aeroelastic wind-tunnel test data 
for the Joined-Wing SensorCraft, a vehicle designed for studying aeroelastic effects using a large number 
of control surfaces and embedded sensors. The method is similar to other real-time frequency response 
methods developed recently 5,6 but mitigates problems associated with spectral leakage and evaluation of the 
Fourier transform at non-integer multiples of the excitation periods. The method also differs from a real-time 
technique that uses a known model structure to estimate parameters, such as stability and control derivatives, 
using the equation-error approach in the frequency domain. 7,8 In addition to computing frequency responses 
from wind-tunnel data, this technique is applicable to identifying aircraft flight dynamics and pilot models, 
aiding in accident investigations, and adjudicating adaptive control laws with simple metrics. 

In section II, the multi-sine input design is summarized, the frequency response and uncertainty estima- 
tion is detailed, and an illustrative example using an F-16 nonlinear simulation is presented. Afterwards in 
section III, the Joined-Wing SensorCraft test article and wind tunnel setup are discussed. In section IV, 
results are presented for bending strain gage measurements due to symmetric control surface deflections. The 
input design, aircraft simulation, and data analysis presented in this work were performed using a MAT- 
LAB toolbox called System IDentihcation Programs for AirCraft (SIDPAC). Developed at NASA Langley 
Research Center, this software is continually improved and is publicly available. 3 
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II. Method 


This section is divided into five parts. The first part summarizes the design of the excitation inputs. The 
second and third parts present the algorithm for estimating the frequency response and its uncertainty from 
measured data. In the fourth part, the method is made recursive for real-time evaluation. The final part 
illustrates the method using an F-16 nonlinear simulation. 

A. Multi-Sine Excitation Design 

Orthogonal, phase-optimized multi-sine perturbation inputs have been developed to excite targeted dynamics 
of multiple- input, multiple-output systems in an efficient manner for system identification. The design of 
these inputs is documented in detail by Morelli 1 ' 9 and briefly summarized here. Each input is the summation 
of sinusoids 

Uj(t) = ^2 a k sin ( ~f~ t + <t>k j (1) 

keKj ' ' 

where Kj is the set of indices for the jth input, a a, is the amplitude, T is the record length, Uk = 2i xk/T is 
the excitation frequency, and 4>k is the phase angle. The record length determines the fundamental frequency 
of the input, of which the excitation frequencies are integer multiples. The excitation frequencies may be 
assigned to the inputs and their sinusoid amplitudes selected arbitrarily. The full set of n w total excitation 
frequency indices K is the union of the disjoint sets Kj. Phase angles are determined using a simplex 
optimization to minimize the peak-to-peak amplitude of each input. 

Superposition of sinusoid inputs minimizes the record length needed for a given amount of excitation 
since all frequencies are excited at the same time, in contrast to traditional frequency sweeps where the 
frequencies are continuously changed. Additionally, the sinusoids are harmonics and therefore orthogonal 
so that all input channels can be excited simultaneously for system identification without correlating the 
responses, which further reduces required data record lengths. Assignment of excitation frequencies to input 
channels and selection of sinusoid amplitudes can be tailored to excite or avoid specific spectral bands and 
to meet signal-to-noise requirements. The optimized phase angles result in small perturbation responses 
from nominal conditions, which are desired for identifying linear models including frequency responses. 
Furthermore, applying all sinusoids simultaneously permits development of the steady-state response, which 
is consistent with the original sine-dwell technique. 

These inputs have led to good modeling results in a variety of unconventional flight conditions such as 
post-stall, high side-slip angle, and hypersonic flight. 7,10 Excitations are added just before software actuator 
limiters on control surface deflection angle and rate to maintain low correlation amidst pilot inputs and 
feedback control. A linear system under this excitation will yield the steady-state responses 

Vi{t) = bksin + (2) 

keK ' ' 

which are sinusoids of the same frequencies but different amplitudes and phase angles. For a system with 
multiple inputs, each output potentially contains dynamic response from all the inputs. 

B. Frequency Response Estimation 

To estimate the frequency response, output measurements are first converted from the time domain into the 
frequency domain. The finite Fourier transform of the output, using the first t seconds of data, is 

y(t,u) = [ y(r )e~ 3WT d.T (3) 

Jo 

where u> is the radian frequency. In applications where the sampling rates are an order of magnitude larger 
than the frequencies of interest, Eq. (3) can be approximated by the discretization 

n — 1 

y(nAt,w) ~ ^ y(mAt)e~ : ’ unAt At, n = 0, 1, . . . ,T/At (4) 

m — 0 

where n is the time index and At is the sampling interval. 
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The Fourier transform of the steady-state response has an analytical solution, found by substituting 
Eq. (2) into Eq. (3), writing the sinusoid in terms of complex exponentials, and integrating by parts. The 
result is 


where 


y(t, u>) = ( fk + jgk) + 6c • (fk - jgk ) 

k£K 

= ^2 + &) ' f k + i{ v k ~ ik) ■ gk 

keK 



(5) 


( 6 ) 


and valid at any time t. Equation (5) is the analytical representation of the steady-state output Fourier 
transform of a linear system subject to multi-sine inputs and is dependent on the excitation frequencies, 
output amplitudes and phase angles, and the record length. If Eq. (5) is evaluated at ui = uik, the formula 
for v k takes the indeterminant form 0/0 and 


lim v k = -jt / 2 (7) 

UJ — 


should be used. 

Equation (4) is used to compute numerical samples of the output Fourier transform, which is represented 
analytically by Eq. (5). The unknown output amplitude and phase angles can be estimated using the 
least-squares framework 

y = X0 (8) 

where y is the vector of rif numerical evaluations of the Fourier transform, 


X = 


(*A+£l) (^2 +&) ••• (lV,+£nJ i j{ v l -£l) 3(^2 - 6 ) ••• j( I V,-£nJ 


(9) 


is an n/ by 2n w matrix of regressors, and 


e = 


/l f'2 ■ ■ ■ fn„ i gi 


lT 


(10) 


is a vector of 2 unknown variables. To keep the least-squares problem well-posed, at least two evaluations 
of the Fourier transform are needed for each excitation frequency; i.e., rif > 2 n u . The least-squares solution 
is 

9 = (11) 

where 3? denotes the real part. At least 2tt/Auj seconds of data are needed so that 3?{XlX} is well-conditioned 
for matrix inversion, where Aw is the smallest increment in excitation frequency. The output amplitude and 
phase angle estimates are extracted from Eq. (11) as 


bk = \/ ft + 9k 

= arctan {gk/ fk)- (12) 

Then given amplitudes and phase angles for the inputs and outputs, the frequency response is computed in 
terms of its relative magnitude and phase angle 

\\G(ju k )\\ = b k /a k 

lG(ju> k ) = i>k~ <t>k ■ (13) 
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Estimation of amplitudes and phases is only required for the output because the input amplitude and 
phase are known from design; however, if actuators introduce significant dynamics or if surface deflection 
measurements are available, the same procedure can be used to estimate amplitudes and phase angles for 
the input. This method mitigates errors from spectral leakage and evaluating the Fourier transform at non- 
integer multiples of the excitation periods because the structure of the steady-state output Fourier transform 
is known. 5 * 


C. Uncertainty Analysis 

Uncertainty in the frequency response is useful for indicating the quality of the estimation. Assuming a 
constant variance of the residual y — X0, the least-squares solution Eq. (11) has the covariance matrix 


S(0) = 


a h 

a fij2 ' 


a h,9l 

a h,92 ■ 

J 1 

a hJl 

a h • 

' °/2,/„„ 

a f 2 ,g i 

a f2,92 ‘ 

J 2 jpn-oj 

Ji 

a f^j2 • 

■ 

a f^ ,91 

a L,92 • 

J riu) 


<,/2 • 

■ a o f 

9 1 ljnu 


2 

a 91,92 ■ 

2 

' <T 9l,9n u , 

°32,/l 

CT 32,/2 ' 

■ a o f 

9 2 ijriuj 

2 

° 92,9l 

CT 92 • 

2 

' tJ 32,3n w 

Ji 

CT L„./2 • 

■ 

a L„,9l 

a L„,92 • 

2 

• <J 9n ul ,9n^ 


= a 2 3?{XTX} 


t-yt- 1 


(14) 


for which the fit error variance estimate is 


<x 2 = 


(y-xe)t(y-xe) 

Uf - 2 n w 


(15) 


and is valid when the problem is over-determined; i.e. , rif > 2n w . 

The uncertainty in the estimation is propagated using the least-squares covariance matrix and first-order 
Taylor series approximations. The output amplitude uncertainty is 


and similarly 





f k ,9k 


— (fk + 9 k) 


C tpk — ( fk+9h ) (dk^fk + fk°g k ~ 2fk9k&f k g k ) 


(16) 

(17) 


for the output phase angle uncertainty. The uncertainties can then be propagated through the the relative 
magnitude and phase angle Eq. (13) as 


a \\G(ju k )\\ 

1 z4 1 

\ k ) 

rr 2 

L G{ju k ) 

= a lk + a lk 



Pki’k ’ 



a o,kbk 


(18) 


The uncertainties of the input amplitudes and phase angles can be neglected because they are either null, 
due to knowledge of the input design, or negligible, as system identification experiments are devised to have 
good measurements of the input. 

Accuracy of error bounds can be improved by obtaining uncertainties on the input and output mea- 
surements, and then propagating the error through the Fourier transforms 11 and applying a weighted least- 
squares solution in place of Eq. (8); however, it was found that this method adds significant computation 
time and was not used. 
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D. Recursive Formulation 


As presented above, the method is used in batch execution after all the experimental data have been collected. 
A reformulation can alternatively be used to update frequency responses in real time as the experiment is 
conducted. The output Fourier transform can be made recursive as 

y(nAt, w) = y(nAt — At, w) + y(nAt)e~^ nAt At. (19) 

For typical record lengths, Fourier transforms may be initiated to zero with negligible losses in accuracy. 
Fourier transforms should be updated with each new data sample, but frequency responses may be estimated 
at any rate desired, provided enough data has been collected to avoid singularities in the least-squares 
estimation. 

E. Illustrative Example 

The method is demonstrated using the flight dynamics of the F-16 fighter aircraft. A nonlinear simulation, 
adapted from Stevens and Lewis, 12 is included in SIDPAC and contains routines for computing linear models 
using central finite differences. The aircraft was trimmed for straight and level flight at 10,000 ft altitude, 
0.24 throttle, 617 ft/s airspeed, and 2.0 deg angle of attack. Numerical linearization resulted in the transfer 
function 

Ag(*) = ~12.10(s + 0.865) 

A S e (s) (s 2 + 2.461s + 7.451) 1 1 

describing the elevator to pitch rate perturbation dynamics with only the longitudinal states and the appli- 
cation of the short period approximation. 13 

The nonlinear simulation was excited with the 20 second nrulti-sine elevator perturbation described 
in Table 1, which has a 2.0 deg peak-to-peak amplitude and frequencies between 0.1 Hz and 2.0 Hz in 
0.1 Hz increments. Figure 1 shows measurements of the elevator, angle of attack, and pitch rate, which 
were corrupted with white Gaussian noise sequences resulting in 20:1, 5:1, and 5:1 signal-to-noise levels, 
respectively. 



Figure 1. Simulated F-16 data for time-invariant modeling 


6 of 17 


American Institute of Aeronautics and Astronautics 



Table 1. Multi-sine F-16 elevator perturbation input parameters (T = 20 s) 


Index 

k 

Amplitude 
afc [deg] 

Frequency 
k/T [Hz] 

Phase Angle 
4>k [rad] 

2 

0.45 

0.1 

1.34 

4 

0.45 

0.2 

0.33 

6 

0.45 

0.3 

4.35 

8 

0.45 

0.4 

3.30 

10 

0.45 

0.5 

6.03 

12 

0.45 

0.6 

2.86 

14 

0.45 

0.7 

4.43 

16 

0.45 

0.8 

6.20 

18 

0.45 

0.9 

2.70 

20 

0.45 

1.0 

2.10 

22 

0.45 

1.1 

6.15 

24 

0.45 

1.2 

0.80 

26 

0.45 

1.3 

4.90 

28 

0.45 

1.4 

2.06 

30 

0.45 

1.5 

5.38 

32 

0.45 

1.6 

5.56 

34 

0.45 

1.7 

0.18 

36 

0.45 

1.8 

0.59 

38 

0.45 

1.9 

0.84 

40 

0.45 

2.0 

4.84 
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The fitting of the input elevator and output pitch rate Fourier transforms at 20 seconds are shown in 
Figures 2(a) and 2(b). The fittings match well, although larger errors in the pitch rate estimation are due 
to increased noise and the presence of the transient response in the Fourier transform. Figures 2(c) and 
2(d) show the magnitude and phase angle estimates and 2cr error bounds at 20 seconds, along with the true 
solution. The frequency response estimates are close to and in statistical agreement with the true frequency 
response. The error bounds are smallest near the short period mode resonance at 2.73 rad/s, where there 
are good excitations and high signal-to-noise ratios. At lower frequencies there are fewer cycles of excitation, 
and at higher frequencies there are lower signal-to-noise ratios. The error bounds are conservative because 
transient solutions are present in the pitch rate response, error in the Fourier transform is not accounted 
taken into account, and the minimal number of regressors are used in the least-squares estimator to save 
computation time. Figures 2(e) and 2(f) show the real-time estimates of the magnitude and phase angle 
along with the final estimates and true frequency response. Although Fourier transforms were recursively 
estimated throughout the maneuver, frequency response estimates were not initiated until 10 seconds to 
ensure that Ji{XfX} was well-conditioned for matrix inversion. Once this condition was met, frequency 
response estimates quickly converged to the true values. 

III. Wind-Tunnel Model and Test Setup 

The Joined-Wing SensorCraft was developed to support the Aerodynamic Efficiency Improvement Pro- 
gram led by the Air Force Research Laboratory. 14 The goals were to demonstrate stable flight under reduced 
static stability, measure stability margins during flight, and alleviate gust loadings on the first two structural 
modes. The model was designed, constructed, and tested by Boeing, the United States Air Force, NextGen 
Aeronautics, Millenium Dynamics, and NASA. 15-20 

The wind-tunnel test article is an aeroelastically-scaled 8% model of the Boeing Joined-Wing SensorCraft 
concept. Bending and torsional stiffness were scaled for the wings and boom using spar design, and the mass 
distribution was scaled using lead weights. The model consists of a fuselage, two aft-swept wings, two 
forward-swept wings, and a tail boom. Figure 3(a) shows a planfornr view schematic indicating locations 
of imbedded instrumentation. There are 13 independent control surfaces which are driven by hydraulic 
actuators and measured using rotational variable differential transformers. Embedded in the wings are 14 
accelerometers, 15 bending strain gauges, 5 leading-edge stagnation point sensors, and 5 rate gyros. 

Several Joined-Wing SensorCraft configurations were tested in the NASA Langley Transonic Dynamics 
Tunnel. The test medium tetrafluoroethane R-134a was used for similitude scaling. The test article was 
mounted on the device shown in Figure 3(b), which afforded optional pitch and plunge degrees of freedom. 
Testing occurred at 0.22 Mach with a dynamic pressure of 51 lb/ft 2 . Data was recorded using a dSPACE 
data aquisition system and decimated to 500 Hz for analysis. Multi-sine excitation perturbations were added 
to trim surface deflection commands using Lab VIEW. 

Highlighted here is test number 614, run number 35, where the pitch and plunge motions were constrained 
and only the aeroelastic deflections were free. Control surfaces were commanded symmetrically about the 
center line. The modeling inputs are the averaged deflection measurements 

u cd = 2^ cRd + $cLd) (21) 

where 6 is the control surface deflection, c £ {FW,AW} denotes forward and aft wing, {R,L} denote right 
and left wing, and d £ {0,M,I} denotes outer, middle, and inner surfaces. Six multi-sine inputs were 
used to excite all 12 surfaces simultaneously, each with a record length of 100 seconds and containing 83 
frequencies between 0.01 Hz and 9.95 Hz. 21 The relatively long duration was needed to accomodate the 
dense frequency gridding desired for aeroelastic analysis. The commanded perturbation deflections had root 
mean squared values of 1.99 deg. Excitation frequencies were distributed in an alternating fashion to the 
input channels and with uniform power. An optimization of the phase angles was prohibitive for such a large 
number of inputs and excitation frequencies, so the phase angles were generated at random and iterated until 
low peak-to-peak amplitudes were obtained. Frequency responses from control surface deflections to sensor 
outputs were estimated; however, in the interest of space, only the forward right wing outer strain due to 
control surface deflection frequency response will be presented here; i.e. , 6 inputs and 1 output. A limited 
time slice sample of model inputs and bending strain gauge outputs is shown in Figure 4. Although ideally 
the output responses should be symmetric about the center line, differences indicate that time histories of 
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(a) input Fourier transform fitting at 20 s 
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(c) magnitude and 2cr errors at 20 s 



(d) phase angle 2cr errors at 20 s 



(e) magnitude evolution 


(f) phase angle evolution 


Figure 2. F-16 elevator to pitch rate frequency response estimation 
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(a) planform schematic (b) testing configuration 

Figure 3. Joined-wing SensorCraft schematic and wind tunnel test setup 20 


the outputs should not be averaged together; rather, each frequency response should be estimated and then 
averaged between the left and right wings. 

As both the actuator command and responses were available for all 12 flaps, system identification was 
performed to acquire actuator models. Figure 5(a) shows a limited time slice sample of the command to 
and response of the forward wing, outer control surface on the left side. Quantization effects are evident 
in the command, as well as amplitude and phase distortions in the response. The equation-error parameter 
estimation method in the frequency domain 3 was used to fit a first-order model to the data, as shown 
in Figure 5(b). Coefficients of determination were above 0.99 and estimated parameter standard errors 
were below 0.2%, implying an excellent match. The spread of identified frequency responses is shown in 
Figure 5(c) and indicates the actuators had a conservative —3 dB bandwidth of 20 Hz, which is twice 
the highest excitation frequency. Despite some magnitude and phase angle distortions, the control surface 
deflections have mutual correlation coefficients below 0.001, suggesting that while they were not orthogonal, 
they were nearly so and served well for system identification. 

IV. Results 

The test article was fixed on the mounting with constrained pitch and plunge degrees of freedom so that 
only the aeroelastic degrees of freedom were excited. The 6 multi-sine perturbation excitations were injected 
at the actuators to deflect the twelve control surfaces in a symmetric manner while sensor measurements were 
recorded. After the test was completed, the real-time frequency response algorithm described was applied 
to estimate the dynamics between the symmetric surface deflections and the outer bending strain on the left 
side forward wing. 

Frequency response estimates obtained at the end of the manuever are shown in Figure 6; forward wing 
surface deflection frequency responses are shown on the left column, and aft wing surface deflection frequency 
responses are on the right column. Consistent with previous findings, the lowest two aeroelastic modes are at 
32.7 rad/s and 63.5 rad/s. 21 The forward wing outer surface forcing had the highest signal-to-noise ratio due 
to the high strains recorded from the co-located actuation and thus the most accurate frequency response 
estimate, despite the fact that all control surfaces were moving simultaneously. Moving inboard along the 
wing and to the aft wing, the frequency responses from the surface deflections had lower strains and lower 
signal-to-noise ratios, resulting in biased frequency response estimates and larger uncertainties. 

The evolution of real-time magnitude and phase angle estimates are shown in Figures 7 and 8. Each trace 
shows the evolution of a magnitude or phase angle estimate at a particular excitation frequency. Frequency 
response estimates were only computed after 50 seconds because that time is needed to ensure that the 
columns of the regressor matrix are linearly independent and the least-squares problem is well-posed. The 
360 deg jumps in the phase angle estimates are an imposed constraint that the lowest phase angle remain 
within [0, 360) deg because this is a causal system. The magnitude estimates of the forward wing outer 
surface deflection to strain converged immediately after 50 seconds, with the exception of estimates around 
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(a) model inputs 

Figure 4. 


(b) bending strain gauge outputs 


Time slice sample of measurements 
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(a) left front wing outer control surface deflection 
command and response 


(b) fit of model to response velocity data using 
equation-error in the frequency domain 




(c) spread of identified actuator frequency responses 

Figure 5. Identified control surface actuator models from command and response measurements 
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the transmission zero at 48 rad/s. Other magnitude frequency response estimates converged in the higher 
frequency range where the signal-to-noise ratios were better. Similarly, the phase angle for the strain due to 
deflection of the forward wing outer surface converged quickly after estimation began, and similarly for the 
other forward wing surface deflections. The phase angle estimates from the aft wing surface deflections had 
lower signal-to-noise ratios and began to converge when the experiment ended at 100 s. 

The forward wing outer surface deflection to outer forward wing strain frequency response was estimated 
well by this method. Frequency responses from the other forward wing surface deflections were estimated 
well only at higher frequencies where the signal-to-noise ratio was higher. Although the structural modes 
were evident in the phase angles, frequency responses from the aft wing surfaces were not estimated well 
because of low signal-to-noise. If this method were applied in real-time as the test was performed, these 
results would be immediately apparent. If only the co- located frequency response was desired, these results 
would be sufficient. If the multi-axis frequency responses were also desired, the test engineers could either 
continue running the tests to accumulate more data, or the amplitudes of the multi-sines could be increased 
for the inner and aft wing surface deflections. It is possible that with this real-time analysis, the results 
could be analyzed as the experiment is conducted to expedite changes in the experiment immediately, rather 
than conducting another test at a later time. This method would save time and money in wind tunnel tests 
by ensuring that results were of desired accuracy. 

V. Conclusion 

In this paper, a new method was presented for estimating frequnecy responses from data in real time 
as an experiment is being conducted. Orthogonal phase-optimized multi-sine inputs were used to compress 
test durations and excite all control inputs simultaneously and uniquely. A recursive Fourier transform was 
used to update the spectral content of each output in a computationally-efficient manner. A least-squares 
method was used to extract amplitude and phase information from these Fourier transforms, which were 
then used to compute the frequency responses and their uncertainties. The method was demonstrated first 
using the short period dynamics of an F-16 nonlinear simulation, where good modeling results were obtained 
using a 20 second maneuver. The method was then demonstrated using Joined-Wing SensorCraft aeroelastic 
wind tunnel data, where frequency responses between 6 symmetric surface deflection inputs and a bending 
strain gauge output were estimated from a single 100 second data record. Good results were obtained for 
the co-located pair, but low signal-to-noise ratios led to poorer results from control surface deflections on 
the aft wing and closer to the aircraft center line. 

A new method was presented for computing frequency responses directly from the data, which bypasses 
many of the engineering judgements needed in conventional frequency response estimation. A measurement 
of the frequency response uncertainty is also produced, which describes the quality of the estimate. Fur- 
thermore, this process can be made recursive for estimation while the data streams in real time rather than 
after the experiment has finished. 

This work is important because it gives direct feedback on modeling results as an experiment is being 
conducted. If adequate results have been obtained, the test can be stopped. If adequate results have not yet 
been obtained, longer test times can be run, larger amplitude inputs can be commanded, or sensor placement 
can be reconfigured. This real-time knowledge will help diagnose problems and expedite wind tunnel tests 
to save time and money. 

This method is useful in aeroelastic tests where there are numerous control inputs that need to be 
excited. The multi-sine inputs allow for simultaneous excitation of all input channels, instead of sequential 
excitation, thereby saving time. The method is useful for aeroelastic wind tunnel tests designed to identify 
frequency responses, which are then used to design control systems for flutter suppression and load alleviation. 
Other applications of this method include identifying aircraft or pilot models and stability margins in flight, 
adjudicating adaptive nonlinear control laws with simple metrics, and alerting pilots and reconfigurable 
control algorithms to changes in aircraft dynamics due to structural damage, loss-of-control situations, and 
ice formation. 
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Figure 6. Final frequency response estimates and 2<r error bounds at 100 s 
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Figure 8. Estimated phase angle evolution 
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